home *** CD-ROM | disk | FTP | other *** search
/ The Very Best of Atari Inside / The Very Best of Atari Inside 1.iso / mint / mntlib25 / _mulsf3.cpp < prev    next >
Text File  |  1992-12-12  |  4KB  |  155 lines

  1. | mjr: not needed on the TT
  2.  
  3. #ifndef    __M68881__
  4.  
  5.     .text
  6.     .even
  7.     .globl    __mulsf3, ___mulsf3
  8.  
  9. __mulsf3:
  10. ___mulsf3:
  11.  
  12. # ifdef    sfp004
  13.  
  14. | single precision floating point stuff for Atari-gcc using the SFP004
  15. | developed with gas
  16. |
  17. | single floating point multiplication routine
  18. |
  19. | M. Ritzert (mjr at dfg.dbp.de)
  20. |
  21. | 7. Juli 1989
  22. |
  23. | revision 1.1: adapted to the gcc lib patchlevel 58
  24. | 4.10.90
  25.  
  26. comm =     -6
  27. resp =    -16
  28. zahl =      0
  29.  
  30.     lea    0xfffffa50:w,a0
  31.     movew    #0x4400,a0@(comm)    | load first argument to fp0
  32.     cmpiw    #0x8900,a0@(resp)    | check
  33.     movel    a7@(4),a0@
  34.     movew    #0x4427,a0@(comm)
  35.     .long    0x0c688900, 0xfff067f8
  36.     movel    a7@(8),a0@
  37.     movew    #0x6400,a0@(comm)    | result to d0
  38.     .long    0x0c688900, 0xfff067f8
  39.     movel    a0@,d0
  40.     rts
  41.  
  42. # else    sfp004
  43.  
  44. | single floating point multiplication routine
  45. |
  46. | written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet).
  47. | Based on a 80x86 floating point packet from comp.os.minix, written by P.Housel
  48. |
  49. |
  50. | Revision 1.2, kub 01-90 :
  51. | added support for denormalized numbers
  52. |
  53. | Revision 1.1, kub 12-89 :
  54. | Created single float version for 68000. Code could be speed up by having
  55. | the accumulator in the 68000 register set ...
  56. |
  57. | Revision 1.0:
  58. | original 8088 code from P.S.Housel for double floats
  59.  
  60. BIAS4    =    0x7F-1
  61.  
  62.     lea    sp@(4),a0
  63.     moveml    d2-d5,sp@-
  64.     moveml    a0@,d4/d5    | d4 = v, d5 = u
  65.     subw    #8,sp        | multiplication accumulator
  66.  
  67.     movel    d5,d0        | d0 = u.exp
  68.     swap    d0
  69.     movew    d0,d2        | d2 = u.sign
  70.     lsrw    #7,d0
  71.     andw    #0xff,d0    | kill sign bit
  72.  
  73.     movel    d4,d1        | d1 = v.exp
  74.     swap    d1
  75.     eorw    d1,d2        | d2 = u.sign ^ v.sign (in bit 31)
  76.     lsrw    #7,d1
  77.     andw    #0xff,d1    | kill sign bit
  78.  
  79.     andl    #0x7fffff,d5    | remove exponent from u.mantissa
  80.     tstw    d0        | check for zero exponent - no leading "1"
  81.     beq    0f
  82.     orl    #0x800000,d5    | restore implied leading "1"
  83.     bra    1f
  84. 0:    addw    #1,d0        | "normalize" exponent
  85. 1:    tstl    d5
  86.     beq    retz        | multiplying zero
  87.  
  88.     andl    #0x7fffff,d4    | remove exponent from v.mantissa
  89.     tstw    d1        | check for zero exponent - no leading "1"
  90.     beq    0f
  91.     orl    #0x800000,d4    | restore implied leading "1"
  92.     bra    1f
  93. 0:    addw    #1,d1        | "normalize" exponent
  94. 1:    tstl    d4
  95.     beq    retz        | multiply by zero
  96.  
  97.     addw    d1,d0        | add exponents,
  98.     subw    #BIAS4+16-8,d0    | remove excess bias, acnt for repositioning
  99.  
  100.     clrl    sp@        | initialize 64-bit product to zero
  101.     clrl    sp@(4)
  102.  
  103. | see Knuth, Seminumerical Algorithms, section 4.3. algorithm M
  104.  
  105.     movew    d4,d3
  106.     mulu    d5,d3        | mulitply with bigit from multiplier
  107.     movel    d3,sp@(4)    | store into result
  108.  
  109.     movel    d4,d3
  110.     swap    d3
  111.     mulu    d5,d3
  112.     addl    d3,sp@(2)    | add to result
  113.  
  114.     swap    d5        | [TOP 8 BITS SHOULD BE ZERO !]
  115.  
  116.     movew    d4,d3
  117.     mulu    d5,d3        | mulitply with bigit from multiplier
  118.     addl    d3,sp@(2)    | store into result (no carry can occur here)
  119.  
  120.     movel    d4,d3
  121.     swap    d3
  122.     mulu    d5,d3
  123.     addl    d3,sp@        | add to result
  124.                 | [TOP 16 BITS SHOULD BE ZERO !]
  125.     moveml    sp@(2),d4-d5    | get the 48 valid mantissa bits
  126.     clrw    d5        | (pad to 64)
  127. 2:
  128.     cmpl    #0x0000ffff,d4    | multiply (shift) until
  129.     bhi    3f        |  1 in upper 16 result bits
  130.     cmpw    #9,d0        | give up for denormalized numbers
  131.     ble    3f
  132.     swap    d4        | (we''re getting here only when multiplying
  133.     swap    d5        |  with a denormalized number; there''s an
  134.     movew    d5,d4        |  eventual loss of 4 bits in the rounding
  135.     clrw    d5        |  byte -- what a pity 8-)
  136.     subw    #16,d0        | decrement exponent
  137.     bra    2b
  138. 3:
  139.     movel    d5,d1        | get rounding bits
  140.     roll    #8,d1
  141.     movel    d1,d3        | see if sticky bit should be set
  142.     andl    #0xffffff00,d3
  143.     beq    4f
  144.     orb    #1,d1        | set "sticky bit" if any low-order set
  145. 4:    addw    #8,sp        | remove accumulator from stack
  146.     jmp    norm_sf        | (result in d4)
  147.  
  148. retz:    clrl    d0        | save zero as result
  149.     addw    #8,sp
  150.     moveml    sp@+,d2-d5
  151.     rts            | no normalizing neccessary
  152.  
  153. # endif    sfp004
  154. #endif    __M68881__
  155.